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Abstract 

In this contribution we develop a cut finite element method with boundary value 
correction of the type originally proposed by Bramble, Dupont, and Thomee in [1]. 
The cut finite element method is a fictitious domain method with Nitsche type en¬ 
forcement of Dirichlet conditions together with stabilization of the elements at the 
boundary which is stable and enjoy optimal order approximation properties. A com¬ 
putational difficulty is, however, the geometric computations related to quadrature 
on the cut elements which must be accurate enough to achieve higher order ap¬ 
proximation. With boundary value correction we may use only a piecewise linear 
approximation of the boundary, which is very convenient in a cut finite element 
method, and still obtain optimal order convergence. The boundary value correction 
is a modified Nitsche formulation involving a Taylor expansion in the normal direc¬ 
tion compensating for the approximation of the boundary. Key to the analysis is a 
consistent stabilization term which enables us to prove stability of the method and a 
priori error estimates with explicit dependence on the meshsize and distance between 
the exact and approximate boundary. 


1 Introduction 

We consider a cut finite element method (CutFEM) for a second order elliptic bound¬ 
ary value problem with Dirichlet conditions. In standard fictitious domain CutFEM the 
boundary is represented on a background grid and allowed to cut through the elements in 
an arbitrary fashion. The Dirichlet conditions are enforced weakly using Nitsche’s method 
[2]. We refer to |3], |1], [5], [6], [7], for recent developments of this approach. See also the 
recent overview paper [8] and [9] for implementation issues. 

Cut finite element methods is one way of alleviating the problem of mesh generation 
and allowing for more structured meshes and associated solvers. For this reason, the 
interest for such methods has increased significantly during the last few years; among 
recent contributions we mention the finite cell method of Parvizian, Diister, et ah mm-, 
the least squares stabilized Lagrange multiplier methods of Haslinger and Renard [12] , Tur 
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et al. [13], and Baiges et al. HI; the stabilization of Nitsche’s method by Godina and 
Baiges [T5|; the local projection stabilization of mnltipliers of Barrenechea and Chonly [16] 
and of Amdonni, Moakher, and Renard im. 

In this contribntion we develop a version of CntFEM based on the idea of bonndary 
valne correction originally proposed for standard hnite element methods on an approximate 
domain in [1] and fnrther developed in [18]. Using the closest point mapping to the exact 
bonndary, or an approximation thereof, the bonndary condition on the exact bonndary 
may be weakly enforced nsing Nitsche’s method on the bonndary of the approximate 
domain. A Taylor expansion is nsed to approximate the valne of the solntion on the 
exact bonndary in terms of the valne and normal derivatives at the discrete approximate 
bonndary. Key to the stability of the method is a consistent stabilization term that, also in 
the case of arbitrary cnt elements at the bonndary, provide control of the variation of the 
fnnction in the vicinity of the bonndary. More precisely, the stabilization ensnres that the 
inverse ineqnality necessary to prove coercivity holds and that the resnlting linear system 
of eqnations has the optimal condition nnmber where h is the mesh parameter, 

independent of the position of the bonndary on the backgronnd grid. 

We prove optimal order a priori error estimates, in the energy and norms, in terms of 
the error in the bonndary approximation and the meshsize. Of particnlar practical impor¬ 
tance is the fact that we may nse a piecewise linear approximation of the bonndary, which 
is very convenient from a compntational point of view since the geometric compntations 
are simple in this case and a piecewise linear distance fnnction may be nsed to constrnct 
the discrete domain. We obtain optimal order convergence for higher order polynomial 
approximation of the solntion if the Taylor expansion has snfficiently high order. In par¬ 
ticnlar, for second and third order polynomials we obtain optimal order error estimates in 
the energy and norms with only one term in the Taylor expansion. Note that withont 
bonndary correction one typically reqnires accnracy in the L°° norm for the ap¬ 

proximation of the domain which leads to signihcantly more involved compntations on the 
cnt elements for higher order elements, see [7]. However, also in the case of no bonndary 
valne correction onr analysis in fact provides optimal order error estimates if the approxi¬ 
mation of the bonndary is accnrate enongh and thns we obtain an analysis for the standard 
cnt hnite element method with approximate bonndary. Finally, we also prove estimates 
for the error both on the discrete domain and on the exact domain. The discrete solntion 
on the exact domain is directly dehned by the method since we may inclnde all elements 
that intersect the nnion of the discrete and exact domains in the active mesh. Even thongh 
some active elements may not intersect the discrete domain the resnlting method is stable 
dne to the stabilization term and no anxiliary extension of the discrete solntion ontside of 
the discrete domain is necessary. We present nnmerical resnlts illnstrating onr theoretical 
hndings. 

The ontline of the paper is as follows: In Section 2 we formnlate the model problem and 
onr method, in Section 3 we present onr theoretical analysis, and in Section 4 we present 
the nnmerical resnlts. 
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2 Model Problem and Method 


2.1 The Domain 

Let be a domain in with smooth boundary dQ and exterior unit normal n. We let 
p be the signed distance function, negative on the inside and positive on the outside, to 
dVl and we let Us{dVl) be the tubular neighborhood {x G : |p(a;)| < 5} of dVl. Then 
there is a constant (5o > 0 such that the closest point mapping p[x) : f/^p (911) —)■ dVL is well 
dehned and we have the identity p{x) = x — p{x)n{p{x)). We assume that 5o is chosen 
small enough that p{x) is a bijection. See [IH], Section 14.6 for further details on distance 
functions. 

2.2 The Model Problem 

We consider the problem: hnd u : H —)■ M such that 

— Am = / in n (2-1) 

u = g on 911 (2.2) 

where / G and g G iL^/^(9H) are given data. It follows from the Lax-Milgram 

Lemma that there exists a unique solution to this problem and we also have the elliptic 
regularity estimate 

||w||r/«+2(o) ^ ||/||m'*(o), s > —1 (2.3) 

Here and below we use the notation < to denote less or equal up to a constant. 

2.3 The Mesh, Discrete Domains, and Finite Element Spaces 

• Let Hq be a convex polygonal domain such that C M'^, where Us{fi) = 

Us{dfl) U n. Let h G (0, ho], be a family of quasiuniform partitions, with mesh 
parameter h, of Hq into shape regular triangles or tetrahedra K. We refer to lCh,o as 
the background mesh. 

• Let H/i, h G (0, ho], be a family of polygonal domains approximating H. To each VLh 

we associate the functions Vh '■ dflh —t \i^h\ = 1, and Qh '■ dQh —t such that if 

Ph{x,g) := X + gi'hi.x) then Ph{x, Qh{.x)) G 9H for all x G dVLh- We will also assume 
that p/i(x,?) G for all x G dVLh and all between 0 and Qh{,x). For conciseness 

we will drop the second argument of ph below whenever it takes the value gh(x). We 
assume that the following assumptions are satished 

h/i := i|^'hi|L“(ao,,) = o(h), h G (0, ho] (2.4) 

and 

\Wh- nop\\Loo(^ 9 ^^) = o{l), hG(0,ho] (2.5) 
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where o(-) denotes the little ordo. We also assume that ho is small enough to guar¬ 
antee that 


hG(o,ho] 


( 2 . 6 ) 


and that there exists M > 0 such for any y G the equation, hnd x G dVLh 

and |<^| < 5h such that 


has a solution set Vh with 


= y 


(2.7) 


caidiVh) < M 


( 2 . 8 ) 


uniformly in h. The rationale of this assumption is to ensure that even if is not a 
bijection its image can not degenerate for vanishing h. 


Choice of During computation, typically the quantities that are easily acces¬ 
sible on dflh are Uh and p. The two choices that are natural for are therefore 

■= nh, Qh '■= with solution to p{j)h{x,c,)) = 0 or := n o p and Qh '■= P- 
Both cases requires the solution of nonlinear equations. The computation of Qh using 
Newton’s method in the first case is substantially less costly than that of n op, since 
the first quantity is a scalar and the initial guess p is more accurate. 


Observe that if := n o p then the mapping ph coincides with p{x). It is therefore 
a bijection and all the above assumptions hold by the properties of the closest point 
mapping. This bijection property does not hold in the general case. However, we 
assume that the equation p{ph{x,‘;)) = 0 has at least one solution for every x G 
dQh and ph may then be identified with the solution of smallest magnitude. As 
an example consider the practically important case where dQh is defined by the zero 
level set of a piecewise linear nodal interpolant of the distance function and we choose 
■= ■n/i, with Hh denoting the normal of dflh- That the associated Ph exists for all 
X G dflh follows immediately from the implicit function theorem: the equation in 
p{x -(- o p) =0 has a solution since p is a bijection and then so does p{x -|- = 0 

since Vp ■ Uh > 0 for h small enough. The assumption (2.8) must clearly hold in 
this case, since if it does not then also p must have a critical point in Us^idfl) (since 
Ph ^ P and the number of solutions is bounded below as h —)■ 0), but this contradicts 
the fact that p is a bijection. Moreover we have the estimates 


Sh < Z!^ Ill'll -no p||l~(8!!») < h 


(2.9) 


• Given a subset a; of Dq) let fChipj) be the submesh defined by 

K:h{ix) = {K elCo,h-.Knuj^^ 
i.e., the submesh consisting of elements that intersect cJ, and let 

= h)KeK.h(u)K 


( 2 . 10 ) 


( 2 . 11 ) 
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be the union of all elements in }Ch{u})- Below the L^-norm of discrete functions 
frequently should be interpreted as the broken norm. For example for norms over N'h 
we have 

IK’IIFh- E (2-12) 

KeJChH 

• Let the active mesh K-h be dehned by 

}Ch:=iCh{nunh) (2.13) 

i.e., the submesh consisting of elements that intersect Qh U ff, and let 

A/h := U ri/i) (2-14) 


be the union of all elements in JCh- 

• Let Vo,h be the space of piecewise continuous polynomials of order p dehned on ICo,h 
and let the hnite element space be dehned by 

14 := Vo^hWh ( 2 - 15 ) 


2.4 Extensions 

There is an extension operator E : Lr^(t4o(12)) such that 

< l4llrfqo), s>0 (2.16) 

see [2nj. For brevity we shall use the notation v for the extended function as well, i.e., 
V = Ev on UsQ^fl). 

2.5 The Method 

Derivation. Let / = Ef and u = Eu be the extensions of / and u from to Usq{VL). 
For u G 14 we have using Green’s formula 

(/, 'y)o;, = (/ + Am, - (Am, m)q^ (2.17) 

= (/ + Am, + (Vm, Vm)o^ - {uh ■ Vm, v)du^ (2.18) 

where we used the fact / + Am = 0 on G, while on \ G we have / + Am = Ef — AEu, 
which is not in general equal to zero. Now the boundary condition m = on dfl may be 
enforced weakly as follows 

(/, v)n^ = (/ + Am, v)n^ + (Vm, Vv)n^ - {rih ■ Vm, v)dn^ (2.19) 

- (uop^ - g op^^Uh - '^v)9n^ + I3h~^(u oph- gop^, v)dn^ 

The positive constant (3 must be chosen large enough to ensure stability, cf. below. 
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Since we do not have access to uopj^ we use a Taylor approximation in the direction Uh 


uopf^{x) ^ Tk{u){x) 






( 2 . 20 ) 


where D l^ is the j:th partial derivative in the direction lyh- Thus it follows that the solution 
to (2.1)-(2(2|) satishes 


(/, = (/ + Am, + (Vm, Vv)n^ - (uh ■ Vu, (2.21) 

- (Tk(u) - goph.Uh- + l3h~^{Tk{u) - gop^, v)dn^ 

- {uoph-Tk{u),nh ■ VM)anft + I3h~^{u o p,^-Tk{u),v)dOh 


for all M G 14 . Rearranging the terms we arrive at 


(Vm, Vm)^^ - {rih ■ Vm, v)oq^ 


- {Tk{u),nh ■ ^v)dnn + l3h~^{Tk{u),v)ann 

+ U + 

-{uopk- Tk{u),nh ■ Vv)anh + I3h~^{u op^- Tk{u),v)aah 
= (/, v)q^ - {go pk, Uh ■ yv)dUf, + I3h-^{g opf,, v)dn^ 


( 2 . 22 ) 


for all V G I 4 . The discrete method is obtained from this formulation by dropping the 
consistency terms of highest order, i.e. those on lines three and four of (2.22). 


Bilinear Forms. We dehne the forms 


ao(M,w) := (Vm, Vm;)^^ (2.23) 

- {uh ■ Vm, w)on^ - {Tk{v),nh ■ 'Vw)an^ 

+ [3h~^{Tk{v),w)anh 

ah{v, w) := ao(v, w) + jh{v, w) (2.24) 

Jh{v,w) := [DI^w])f (2.25) 

FeXh ^=1 

lh{w) := (/, w)n^ - {go Ph, nh ■ '^w)ann + l3h~\g o p^, w)an^ (2.26) 


where 7 ^- is a positive constant. Here we used the notation: 

• tFh is the set of all internal faces to elements K E K.h that intersect the set Q\Qht-)dQh, 
and np is a. hxed unit normal to F G J4- 

• is the partial derivative of order I in the direction of the normal rip to the face 

FeJ^h. 

• Hlf’ = with = lims^o+ v{x^snp)^ is the jump of a discontinuous function 

V across a face F E Fh- 
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The Method. Find: Uh G Vh such that 


ahiuh.v) = Ihiv), 


where ah is dehned in (2.24) and //, 




Vn e Vh 


(2.27) 


Symmetric Formulation in the Case fc = 1. Using one term in the Taylor expansion 
gives the following forms 


ah{v, w) = (Vn, + jh{v, w) (2.28) 

- {uh ■ Vn, w) 9 n^ - {v, Uh ■ 'Vw)on^ 

- {Qh^h ■ 'Vv,nh ■ 'Vw)an^ 

+ /3h~^{Ti{v),w)an^ 

h{w) = (/, w)n^ - {go Ph, rih ■ '^w)an^ + /3h-^{g o w)an^ (2.29) 


We see that only the terms of the third and the fourth lines of (2.28) violate the symmetry 
of the formulation. To make it symmetric we choose Uh := rih, assuming that the discrete 
approximation Qh is such that this is a valid choice and also symmetrize the penalty term 
in the fourth line by replacing w in the right hand slot by Ti{w). A similar perturbation 
is added to the right hand side to keep consistency. The forms of the resulting symmetric 
formulation read 


ah{v, w) = (Vn, V«;)nh + jh{v, w) (2.30) 

- {uh ■ Vv,w)an,, - {v,nh ■ 'Vw)anh 

- {Qh rih ■ Vn, Uh ■ 'Vw)anh 

+ ^h-\T,{v),T,{w))an, 

h{w) = (/, - {go Ph, Uh ■ Vw)aQf, + I3h~^{g o ph, Ti{w))an^, (2.31) 

The analysis presented below covers this important special case. Also observe that if 
more terms are included in the Taylor development the resulting nonsymmetric part of the 
matrix is expected to be small, relative to the symmetric part, and the reduced symmetric 
form is likely to be a good preconditioner. 


3 A Priori Error Estimates 

3.1 The Energy Norm 

Let the energy norm be dehned by 

lll^lll^ = llVt^ll^, + lllt^llll, + hWrih ■ VnllL, + h-^v\\l^^ (3.1) 

where 

IMW^ = jh{v,v) (3.2) 
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3.2 Consistency 

In view of (2.22) we obtain the identity 


ah{u - Uh, v) = {uoph- Tk{u),nh ■ - /3h ^{uopf^- Tk{u),v)gn^ (3.3) 

+ (/ + An, v)n^\n, Wv G 14 

and thus we conclude that 

\ah{u - Uh,v)\ < \\uopf, - Tk{u)\\on^(j\nh ■ + h~^/3\\v\\on^'^ (3.4) 

+ \\f + 

< h~^^‘^\\uopi^-Tk{u)\\gnJ\\v\\\h (3.5) 

+ II/ +An|4^\n||'y||o^\n, Vn e 14 


Estimate of the Error in the Taylor Approximation. The Taylor polynomial 
Tk{u){x) provides an approximation of u opj^{x) and we have the error estimate 


\vopj^{x) - Tk{v){x)\ 


< 




Dlpv{x(s)){e^(x) - s'fds 




(3.6) 

(3.7) 

(3.8) 


where 4 is the line segment between x and Ph{x). Combining (3.4) and (3.8) and recalling 


the assumption (2.8) we arrive at the estimate 


\v o p - Tk{v)\\l^^ 


< 

f 1 

\Di:^v\ 

|2 

\lx 

\gh{x)\‘^’'^^dx 

(3.9) 


JaUh 




< 

rsj 

/ 1 

\di:m 

|2 

\ls, 

\gh{x)\‘^’'^^dx 

(3.10) 





< , 

rsj 


\\d’^+\\ 

\\2 

m 


(3.11) 


Here we handled the possible overlap of the contributions from different polygonal sides of 


dilh by using the fact that by assumption (2.8) such an overlap must have a hnite num¬ 


ber of contributions uniformly in h and by dropping the directional derivative, effectively 
including the derivatives of order fc + 1 in all directions. 

With slightly stronger control of the regularity we obtain the estimate 


\v o p - Tk{v)\\gn, 


<5^+^ sup 
0<t<So 


(3.12) 


where dflt = {x G : p{x) = t} is the levelset with distance t to the boundary dfl. 









Estimate of the Residual on \ Suppose that 

f + Au € H‘(Us,(,a)) 


(3.13) 


which, in view of (2.3) and (2.16), holds if / G Using (3.13) and the fact that 

/ + Am = 0 in fl, we obtain the estimate 


ll/ + A'“lb.\n;S<l|i>»(/ + A'“)llnAn:£C'" l|£'U/ + A'«)||an, 

0<t<So 

where we used the fact that fl/i \ U C Us{dfl), where 6 ^ 6h- 


(3.14) 


Estimates of the Consistency Error. Combining (3.12), (3.14), and (3.16), we obtain 
the estimate 

\ah{u-Uh,v)\ <6^+^ sup \\D’'+^u\\L^ant)(\\nh-'^v\\on^ + h-^/3\\v\\onA (3.15) 

0<t<5o ^ ^ 

+ sup ||U>^(/ +AM)||anJ|'y||n;,\o, Vn e 14 

0<t<So 

This estimate will be used when we derive an estimate of the error while for the energy 
error estimate we continue the estimation using the bound 


which leads to 

\ah{u-Uh,v)\ < sup 42(9^0 

^ 0<t<5o 

+ 54^ sup ||T>^(/ +AM)| 4 nt)|| 4 IIU, ^veVh 

0<t<So ' 


(3.16) 

(3.17) 


Remark 3.1 VTe may upper bound the right hand sides further using global trace inegual- 
ities leading to 

sup < \\u\\Hk+ 2 ^Q) < ||/||m'=(o) (3.18) 


and 


0<t<5o 


sup \\Dhif + Am)|4qj < ||/||/f!+i(o) + IIAM|4i+i(f7) < ||/||m*+i(o) 
0<t<5o 


(3.19) 


3.3 Inverse Inequality 


Using the additional stability provided by the stabilization term we have the following 
inverse inequalities 

l|Vn||^^<||VM|4^ +114111,\, VneU, (3.20) 

and 

l|u|lv, 4 Iklln^+ h2||4||||,, VneU* (3.21) 

See [6] for a proof. 
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3.4 Coercivity and Continuity 

We have coercivity 

\M\\l^ahiv,v), Vhe(0,ho] 
if ho small enough and /? large enough, and continuity 


(3.22) 


(^h{v,w) < (^\\\v\\\h +h |||m;|||/i, Vn e Id + 14, tc e 14 (3.23) 


where 


1 / = //^+^/2(A4) n n 


(3.24) 


is the space on which the functional V 3 v ^ ah{v,w) G M, for a hxed w & Vh and 
hxed h G (0, ho] is bounded. The continuity estimate ( |3.23 ) follows directly from the 
Cauchy-Schwarz inequality and we next verify the coercivity estimate (3.22). 


Verification of (3.22). Using the notation 

Ti,k{v) = Tk{v) - V 

we obtain 

ah{v, v) = (Vn, Vv)n^ -h jh{v, v) - 2{nh ■ Vn, n)an^ (3h~^{v, n)ao^ 
/3h~^(Ti^k{v),v)on^ - (Ti^k{v),nh ■ Vn)an^ 

^ liVnIln;, + lll'f^llli - h^^^Wnh ■ Vv\\on^h-^^'^\\v\\on^ + ^h-^\ 


IdUh 


(3.25) 


(3.26) 


(3.27) 


- Ph-^/^T,,k{v)\\dn,h-^^^\\v\\9n, - h-^/^\\T,,k{v)\\an,h^^^nk ■ 
Now we have the inverse bounds 


h^/^Wuh ■ Vnllan, < ||Vt;||Ar,(9o,) 


(3.28) 


/^-'/^||Ti,.(n)||ao.<h-iTi,,(n)|U,(ao, 


< 

rs_/ 


si 




0=1 


^o(h) 

< 7(^)l|Vn||Ara90;,) 


(3.29) 

(3.30) 

(3.31) 


where 7 (h) —)■ 0 as h —)■ 0. Using (3.20) together with obvious estimates these bounds 


conclude the proof of the coercivity result (3.22) for [3 large enough and h G (0, ho], with 
ho small enough. 


□ 


10 










3.5 Interpolation Estimates 

Let 

TTh ■ 3 M I—)■ 'KsZ,hEu G Vh 


(3.32) 


where E is the extension operator introduced in Section 2.4, and T^sz.h is the Scott-Zhang 
interpolation operator. The following error estimate for the Scott-Zhang interpolant is well 
known 


\\u-'Ksz,hu\\Hm^K)<h" ^<m<s<p + l, K e tCh (3.33) 

Using the properties of the extension operator we then immediately deduce this interpola¬ 
tion error estimate for (3.32) 


\u - 7ihu\\\h + h - 7ihu)\\anh ~ hP\\u\\HP+^n) 


(3.34) 


Verification of (3.34). The first term in (|3.34) is estimated using the trace inequality 


1(90, 


\\v\\K + h\\Vv\\K, 


K elCh 


(3.35) 


see [22], followed by the interpolation estimate (3.33) and stability of the extension operator 
(2.16). Again using the trace inequality (3.35) the second term in (3.34) can be estimated 
as follows 


h - vr/,n)||ao„ < h ^\\Ti^k{u - 7r/jM)||Ar„(ao„) + ||VTi,fc(n - (3.36) 

< hP\\u\\HP+.^n) (3.37) 

where finally we used the fact that Sh ^ h and the estimate 


k 


h”^-i||V’"Ti,fc(n-7r,n)||x 

3=^ 

(3.38) 


k 

7 = 1 

(3.39) 


J ^ 

< m| J7P+i(At(iC)) 

(3.40) 

for m = 0,1 and K e K.h{dVLh)- 


□ 


3.6 Error Estimates 

Theorem 3.1 If Sh = o{h), then the following estimate holds 

\\\u-Uh\\\h<hP\\u\\HP+^(Q.) + h~^/‘^5l"^^ sup (3.41) 

0<t<So 

+ 6^+^ sup \\Dl{f + Au)\\L 2 (^g^^) 

-So<t<0 
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p 

hP 

k 


1 

sr 

1 

h^ 

0 

hl.T> 

0 

h'^ 

2 

h^ 

1 


1 

h^ 

3 

h^ 

2 


2 


4 

h^ 

3 


3 



Table 1: The order of the terms in the energy error estimate under the assumption 5h^h?- 
We conclude that we obtain optimal order of convergence for p = 2,3, with one term, k = 1, 
in the Taylor expansion and for p = 4, 5, with two terms, k = 2. 


Proof. We hrst note that adding and subtracting an interpolant and using the triangle 
inequality and the interpolation estimate (3.34), we obtain 


\\\u - Uh\\\h < |||m - Tihu\\\h + WlnhU -Uh\\\h (3.42) 

< hP\\u\\HP+^{n) + WlTThU - Uh\\\h (3.43) 

For the second term on the right hand side we have the estimates 

11 K/iM -Uh\\\l< tthillhU - Uh, TThU - Uh) (3.44) 

= tthiTThU - U, TThU - Uh) + ahiu - Uh, TThU - Uh) (3.45) 

< {\\\TThU-u\\\h + h~^/‘^\\Ti^k{T^hU-u)\\3^^}j\\\TThU-Uh\\\h (3.46) 

+ h~^/‘^\\uoph - Tk{u)\\3<^^\\\TThU - Uh\\\h 

+ 11 / + 

< h^\\u\\HP+^(n)\\\T^hU - Uh\\\h (3.47) 

+ sup \\D^^^u\\Lii^Qnt)\\\TThU - Uh\\\h 

0<t<(5o 

+ sup + /\u)\\L 2 {dnt)\\\T^hU- Uh\\\h 

-(So<t<0 


where we used coercivity (3.22), added and subtracted the exact solution u, estimated the 
hrst term using continuity (3.23) followed by the interpolation estimate (3.34) and the sec¬ 
ond using the consistency estimate (3.6). Combining estimates (3.43) and (3.47) concludes 
the proof. □ 


Theorem 3.2 


If5h< 


h?, then the following estimate holds 

l|e|b;. < hP+^\\u\\HP+^n) 

+ 5^+^ sup \\D^+\\\L2(^ant) 
0<t<So 

+ sup \\Dlif + Au)\\on, 

0<t<So 


(3.48) 
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p 

h-P+i 

k 


1 


1 


0 


0 


2 


1 


1 


3 

h* 

2 


2 


4 


3 


3 



Table 2: The order of the terms in the energy error estimate under the assumption that 
6h ^ We conclude that we obtain optimal order of convergence for p = 2,3, with one 
term, k = 1, in the Taylor expansion and for p = 4,5, with two terms, k = 2. 


Proof. Let 0 G Hq{Q) be the solution to the dual problem 

a{v,4)) = {v,il)), V & Hl{Vt) (3.49) 

where ip = u — Uh on klh and = 0 on 12 \ ff/i, and extend p using the extension operator 
to Us^lkl). Then we have the stability estimate 

11011^2(0) < Il'i/’llonOh (3.50) 

We obtain the following representation formula for the error 

llello;, = {e,P + A0)o^ - (e, (3.51) 

= {e,p + A0)Q^\n + (Ve, Vp)n^ - (e, Uh ■ V0)an^ (3.52) 

= {e,p + Ap)n^\n + ao(e, p) + bh{e, p) (3.53) 

= / + // + /// (3.54) 

where 

III = (Tfc(e) -e,nh- Vp)dnh - /3h~^{Tk{e),p)dn^ + {uh ■ Ve,p)dQ^ (3.55) 

= (Ti,fc(e), n ■ Vp)dn^ - /3h~^{e, p)dn^ (3.56) 

- /3h~^{Ti^k{e), p)dnh + i'^^h ■ Ve, p)an^ 


Term I. We have 

\I\ = \ie,P + Ap)n,\n\ (3.57) 

^ ||e||o;,\Q||^ + A0||Q^\n (3.58) 

< ({Sl + h 6 h)\\\e\\\l^ (^||e||o;^\o + ||0||//2(Q)j 

< {h~‘^6h + /i"^4)^'^^^/i|||e|||h||e||n;, 

'V 
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(3.59) 

(3.60) 

(3.61) 




Here we used the estimate 


+ veH^Qh) (3.62) 

with V = e, the dehnition of the energy norm to conclude that h“^||e||gQ^ < ll|c|||^, the 
stability (2.16) of the extension operator, the stability (3.50) of the dual problem and the 
assumption that 5h^h?- 


Term II. Adding and subtracting an interpolant we obtain 


n\ = |aft(e, 0 - 7rfe0) + a/i(e, | (3.63) 

< |||e|||ft|||0 - + |a/,(e,7r/j0)| (3.64) 

^ ^ll|e||U||0||r/2(o) + \ah{e,'Kh4>)\ (3.65) 

^ ^ll|e||UI|e||o^ + \ah{e,Tih(t>)\ (3.66) 


To estimate the second term on the right hand side we employ (3.15), with v = 7r/j0, 


|afe(e,7r/j0)| < sup \\D^""^u\\L2i^dnt)(\Wh-^T^h(t)\\ann +h ^(5\\nh(t)\\dnA (3.67) 

0<i<5o ^ ' 

+ sup \\Dl{f + Au)\\ant\Wh(j)\K\n 

0<t<(5o 


Here we have the estimates 


\\nh ■ V7r0||an^ + h ^\\7i(l)\\aQ^ < \\nh ■ V(7r0 - 0)||ao^ + h ^||7r0 - (l)\\an^ (3.68) 

+ \\nh ■ V0||ao^ + h~^\\(j)\\anh 

< h-^/^\\\7i(j) - (j)\\\h (3.69) 

+ HWH^Uh) + h~^^h‘^MH^{Us^(an)) 

^ h^^‘^\\4>\\H^(n) + Il0l|rr2(n^) + h ^Sl^‘^\\(j)\\H^Us^^an)) (3.70) 

< (h^/^ + 1 + h ^Sh‘^)\\e\\n^ (3.71) 

< l|e|b. (3.72) 


and 


||7^/^0||nA^ < + H\K\^ (3-73) 

< h‘^\\(j)\\H2{n) + Sh\\'^(j>\\usf^{dn) (3.74) 

< {h^ + SMK (3.75) 

< WK (3.76) 

where, in both estimates, we used the assumption 5h^h?, as well as the following bounds 

l|.^IUn,. < sfWn ■ (3.77) 
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||0lbh\n ^ \\4>\\us.{an) ^ Sh\\n-V^UsAan) 


(3.78) 


see the Appendix for the proof of these estimates. Combining estimates (3.67), (3.68), and 


(3.73), we arrive at 


|afe(e, 7 r/j 0 )| < (5^+^ snp \\D^^\\\L 2 (on,) 

^ 0<t<So 

+ s‘A" sup iia‘.(/ + a«)|| aut') lleiloft 

0<t<So ^ 


(3.79) 


which together with (3.66) gives 


I//I < 


e\\\h + S^'^^ snp ||A>^+^n||L 2 (aop 
0<t<5o 

snp ||A)(,(/ +AM)||aQ,)||e||n, 

0<t<5o ' 


Term III. Using the Cauchy-Schwarz ineqnality we get 

\IU\ = \hh{eM 

~ \\Ti,k{e)\\aQ.h\\nh ■ V(/)||an;, + /3h"^||e||anhll0ll(9Oh 
+ /?h-i||Ti,fe(e)|| + W'^h ■ Ve||anhll0||an,i 

< \\Ti,k{e)\\dnn (^h"^||0||ao;, + \\nh ■ V0||an;, 

+ |||e|||hh ^'^^ll^llaOh 

^ (^||7'i,fc(e)||aQ;, + h~^^‘^Sh\\\e\\\h^ ||e||oft 


(3.80) 


(3.81) 

(3.82) 

(3.83) 

(3.84) 

(3.85) 


where we used (|3.77|) and (|3.78|) followed by the stability estimate for the dual problem 

(3.86) 


(3.50), and at last the estimate 

||7'i,fe(e)||ao^ < h^~^^\\u\\HP+An) + {h~^^‘^6h)h\\\e\\\h 


Verification of (3.86). We have 


l|r.,»(e)|| 8 u,.< 5 ;ii||£>J,e|| 8 n. 

i=i 


(3.87) 


and for each of the terms ||Zl^^e||ao^, j = 1 ,..., fc, we obtain by adding and subtracting 
an interpolant, using the interpolation estimate (3.33) for the first term and an inverse 
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estimate for the second, the estimates 

^ h~^\\Di^iu - 7ihu)\\]^^^Qn^) + h\\VDi^{u - 

+ h-^\\Di^{T^hU - UhWj^^^en^^ + h\\VDi^{7rhU - UhWj^^^Qn^^ 
< h2^+^"^^'||M||ip+i(Ar,(At^(ao,))) + h^~^mv{7rhu - 


which leads to 

sm\Di -"^ 


Uh^WdUh ~ 


< h-'(S,/h)yh^^’’+»\\up„,„^„^ + A(4/A)"''l|Ve||^^,8n^, 

where we used (2.4) and the fact 5h/h? < 1. Thus we have 

k 


l|ri,^WI|an, <5;^il|C4e||an. 

i=i 


<{h-^/^6,)(hP+^\\u\\H.^r^n) + 


eh 


(3.88) 

(3.89) 

(3.90) 

(3.91) 

(3.92) 

(3.93) 


(3.94) 


Conclusion of the Proof. 

I, II, and III, we obtain 


Collecting the bounds (3.61), (3.80), and (3.85), of Terms 




< 

rsj 


e\\\h 


+ h^’'’^||M||r/P+i(n) 

+ (5^+^ sup \\D’''^\\\L 2 (^on,) 
0<t<So 

+ sup \\Dl{f + Au)Won, 


0<t<5o 


which together with the energy norm error estimate (3.41) concludes the proof. 


(3.95) 


□ 


Theorem 3.3 

and 


The following estimates hold 

ll^^lln ^ h^\\u\\HP+^{n) + |||e|||h 

Iloilo ^ h^~^^\\u\\HP+i(n) + ||e||n^ + h|||e|||/i 


(3.96) 

(3.97) 


Proof. Adding and subtracting an interpolant, using the interpolation estimate (3.33), 
and the inverse inequality (|3.20[) or (3.21), we obtain, for m = 0,1, 


II < II V”^(« - 7r,«)|b\o, + II V”^(7r,« - u,)\\n\n, (3.98) 

< ||V™(7r/jM - Uh)\\n^ + h^~"'\\\7rhU - Mh|||j^ (3.99) 

< hP+^-^\\u\\HP+^n) + ||V™e|b, + hi—lllelll,, (3.100) 

< hP^^-^\\u\\Hp+Hn) + ||V™e|b, + h^—|||e|||, (3.101) 
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which concludes the proof. 


□ 


Remark 3.2 If for a given p the lowest values of k and I are chosen so that optimal 
convergence is obtained, it is straightforward to use a trace ineguality, see (3.18) and (3.19), 
to show that 


\\u — + h|||rt — M/i||| < h^~^^{\\f\\HP-i{n) + ||M||rfp+i(o)) 

Therefore the regularities reguired for optimality of the consistency error of the boundary 
approximation are always optimal compared to the polynomial approximation. 

Remark 3.3 We note that we obtain, as a special case, optimal order error estimates for 
the standard cut Nitsche method with approximate domains by assuming k = 0 and 

6h < (3.102) 

for the energy norm estimate and 

5h < (3.103) 

for the if norm estimate. The latter assumption is comparable with the geometric approx¬ 
imation accuracy achieved by standard isoparametric finite elements of order p. 


4 Numerical Examples 

In the numerical examples, we use implicitly defined boundaries by use of zero isolines 
to predefined functions. Two examples have been considered, one with both convex and 
concave boundaries, so that cut elements can have parts outside the actual domain, and 
one example with nonzero boundary conditions where we also compare setting the bound¬ 
ary condition on the exact boundary to setting them on computational boundary. In all 
examples the stabilization parameters were set to 7 ^ = 1 / 10 , /d = 100 . 


4.1 Convex and Concave Boundaries 

In our first exampe we consider a ring-shaped domain. In Fig. [^we show the zero isoline 
of the function cf = (i? — 1/4) (R — 3/4), R = \/x'^ -\- y'^, used to implicitly define the 
domain, and the resulting mesh after removing the cut part. On this ring, we used a load 
corresponding to the exact solution being a square function in R, 

M = 20(3/4-R)(R-1/4) (4.1) 


with zero boundary conditions on the outside as well as inside boundaries. The elements on 
the inside of the ring are partially outside the computational domain; outside the domain 
the load was extended by zero and the exact solution (in the convergence study) by (4.1). 
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We show an elevation of the approximate solntion on one of the meshes in a seqnence in 
Fig- In Fi g. [3] andj^we show the convergence rates obtained using the symmetric method 
(2.30)-(2.31) for and elements (polynomial orders p = 2 and p = 3), respectively. 
We also show the suboptimal convergence rates of the original Nitsche method. Note in 
particular that the optimal rate is attained also for p = 3 even though only the hrst two 
terms in the Taylor series are accounted for. 


4.2 Nonzero Boundary Conditions 

The domain for the second exampe lies inside the ellipse dehned by the zero isoline to 
0 = /{3/AY + |/^/(l/2)^ — 1. In Fig. [^we show the zero isoline of this function and the 

resulting mesh after removing the cut part. On this domain we use the right-hand side 

/ = TT^ cos ( 7 rx/ 2 ) cos {Tiy/ 2 ) 


corresponding to the exact solution u = cos (7ra;/2) cos (7r|//2). This function also dehnes 
the boundary conditions on the cut boundary. An elevation of an approximate solution on 
one of the meshes in a sequence is given in Fig. 

In Fig. we show the observed L 2 convergence with a P^ approximation using four 
different approaches: 


The symmetric method (2.30)-(2.31). 


• The unsymmetric Taylor expansion with two terms. 


• The unsymmetric Taylor expansion with three terms. 

• Prescribing the boundary condition on the cut boundary (using the fact that the 
exact solution is known). 


In all cases the rate of convergence is 4, which is optimal. The error constant is slightly 
better if we prescribe the boundary condition on the cut boundary, which is to be expected 
since this does not introduce any approximations of the boundary condition. The difference 
between the other three methods is negligible. 


Appendix: Verification of Some Estimates 

Estimate ( |3.77D . For each x G dVth we have the representation 

4){x) = 4){p{x)) + / V(j){sx + {1 — s)p{x)) ■ {x — p{x)) ds (4.2) 

Jo 

Using the Cauchy-Schwarz inequality we obtain 

l0(^)P<l0(p(^))P + 4||n-V0||l (4.3) 

< |0(p(x))|2 + 6h\\n ■ (4.4) 
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where Ix is the line segment between x and p{x) and hf^{j){x)) is the line segment between 
the points p{x) ± 5hn{p{x)). 

Integrating (4.4) over we obtain 

^ I / X 11 .^ ' 


II0II 




'dVlh 


'dn^ 


< 


< 


IdQ 


\(l>{y)\dy+ dh\\n-V(j)\\r^dy 


'an 


+ dh\\n ■ (an) 



(4.5) 

dy 

(4.6) 


(4.7) 


where we first changed the domain of integration from dQh to dQ and then from the 
tubular coordinates to the Euclidian coordinates. Next integrating (4.4) over Is^{y), with 
y = p{x) G dfl, we obtain 




Shiv) 


Again using appropriate changes of coordinates we obtain 


\\^\\us^(aQ) ~ 


< 


'an 


MlAy)dy 


'an 


Sh\(l){y)\‘^dy+ I Sl\\n ■ Vm^^^^dy 


'an 




(4.8) 

(4.9) 

(4.10) 

(4.11) 
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Figure 1: Background mesh with the boundary of indicated, and the corresponding 
computational mesh. 


22 










Figure 2: Elevation of the approximate solution on one of the meshes in a sequence. 
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Figure 3: Convergence using elements, symmetric form. 
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Figure 4: Convergence using elements, symmetric form. 






Figure 5: Background mesh with the boundary of indicated, and the corresponding 
computational mesh. 
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